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ABSTRACT 

-  A  large-eddy  simulation  (LES)  of  transition  in  plane  channel  flow  has  been  carried  out. 
The  LES  results  have  been  compared  with  those  of  a  fine  direct  numerical  simulation  (DNS), 
and  with  those  of  a  coarse  DNS  that  uses  the  same  mesh  as  the  LES,  but  no  residual  stress 
model.  While  at  the  early  stages  of  transition  LES  and  coarse  DNS  give  the  same  results, 
the  presence  of  the  residual  stress  model  was  found  to  be  necessary  to  predict  accurately 
mean  velocity  and  Reynolds  stress  profiles  during  the  late  stages  of  transition  (after  the 
second  spike  stage).  The  evolution  of  single  Fourier  modes  is  also  predicted  more  accurately 
by  the  LES  than  by  the  DNS.  As  small  scales  are  generated,  the  dissipative  character  of 
the  residual  stress  starts  to  reproduce  correctly  the  energy  cascade;  as  transition  progresses, 
then,  and  the  flow  approaches  its  fully  developed  turbulent  state,  the  subgrid  scales  tend 
towards  equilibrium  and  the  model  becomes  more  accurate. 
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1.  Introduction 

In  large-eddy  simulations  (LES)  of  the  Navier-Stokes  equations  the  large,  energy-carrying 
scales  of  motion  are  accurately  captured,  while  the  effect  of  the  small  scales  is  modeled. 
The  models  used  to  parameterize  the  effect  of  the  subgrid  scales  (known  as  residual  stress 
models  or  subgrid-scale  stress  models)  are  simpler  and  more  universal  than  the  turbulence 
models  used  for  the  Reynolds- averaged  Navier-Stokes  equations,  given  that  the  subgrid  scales 
tend  to  be  more  isotropic  and  homogeneous  than  the  large  ones,  and  less  affected  by  the 
boundary  conditions.  Since  only  the  large  structures  are  accurately  resolved  by  the  grid,  LES 
requires  less  computational  effort  than  the  direct  numerical  simulation  (DNS)  of  the  Navier- 
Stokes  equations,  in  which  all  scales  of  motion  are  accurately  resolved,  at  the  expense  of 
the  empiricism  introduced  by  the  modeling  assumptions.  Thus,  one  can  perform  large-eddy 
simulations  of  a  given  flow  at  a  fraction  of  the  expense  required  by  the  DNS,  or,  conversely, 
one  can  study  higher  Reynolds  number  flows  by  LES  than  one  can  by  DNS. 

Large-eddy  simulations  have  been  successfully  applied  to  a  variety  of  turbulent,  wall- 
bounded  flows  such  as  plane  channel  flow  [1,  2,  3],  boundary  layers  [4]  and  channel  flow  with 
transpiration  [5],  but  only  recently  efforts  have  been  made  to  study  transition  to  turbulence 
using  LES  [6,  7,  8,  9].  While  earlier  work  [6,  7]  was  characterized  by  the  application  of  well- 
established  residual  stress  models  to  the  simulation  of  laminar-turbulent  transition,  Piomelli 
and  co-workers  [8]  were  the  first  to  use  the  databases  generated  by  the  direct  simulation 
of  transition  in  a  plane  channel  to  study  the  behavior  of  the  residual  stress  tensor  during 
transition.  They  observed  that  during  the  nonlinear  stages,  and  in  particular  during  the 
second-spike  stage,  the  subgrid  scale  dissipation  (f.e.,  the  energy  transfer  from  large  to  small 
scales)  is  significantly  smaller  than  in  turbulent  flow.  They  devised  a  modification  of  the 
Smagorinsky  [10]  model  which  allowed  accurate  prediction  of  tne  early  stages  of  transition 
in  a  flat-plate  boundary  layer.  Subsequently,  in  [9],  they  again  calculated  a  transitional 
boundary  layer  flow  using  a  residual  stress  model  based  on  the  Renormalization  Group 
(RNG)  Theory  [11]  which  predicts  zero  eddy- viscosity  as  long  as  the  magnitude  of  the  strain 
rate  tensor  is  less  than  some  threshold  value. 

The  latter  results  [8,  9],  however,  indicate  that,  up  to  the  early  stages  of  the  laminar 
breakdown,  the  effect  of  the  small  scales  on  the  resolved  ones  is  negligible;  a  coarse  direct 
simulation,  then,  appears  to  be  sufficient  to  give  results  in  satisfactory  agreement  with 
those  of  finely  resolved  direct  simulations  [12],  as  far  as  engineering  statistics  are  concerned 
(Reynolds  stresses  and  mean  velocity  profiles,  for  example).  For  this  reason  the  present 
article  will  focus  on  the  late  stages  of  transition,  and  apply  the  modified  Smagorinsky  model 
[8]  to  the  study  of  transition  in  a  plane  channel.  The  present  calculation  will  extend  through 
the  entire  transition  regime  and  into  turbulent  flow.  Large-eddy  simulation  results  will  be 
compared  with  those  of  a  finely  resolved  direct  simulation  and  with  those  of  a  coarse  direct 
simulation  to  distinguish  modeling  from  resolution  errors,  and  ascertain  whether  or  not  the 
presence  of  a  residual  stress  model  is  required  for  accuracy. 

In  Section  2  the  governing  equations  will  be  shown.  In  Section  3,  the  results  of  the 
large-eddy  simulation  of  transitional  channel  flow  will  be  discussed.  Conclusions  and  recom¬ 
mendations  for  future  work  will  be  made  in  the  last  section. 
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2.  Numerical  formulation 


In  large-eddy  simulations  we  decompose  the  flow  variables  into  a  large  scale  (or  resolved) 
component,  denoted  by  an  overbar,  and  a  subgrid  scale  component.  The  large  scale  compo¬ 
nent  is  defined  by  the  filtering  operation: 


/(xi ,  ^2,  ass)  —  Gi(xj,  2;j)/(xj ,  ^2,  X2^dx-^dx2dx2 


(1) 


where  the  integral  is  extended  over  the  entire  domain  D  and  Gi  is  the  filter  function  in  the 
ith  direction  [in  the  present  work  Xi  or  x  is  the  streamwise  direction,  xj  or  y  is  the  direction 
normal  to  the  walls,  which  are  located  &i  y  =  ±1,  and  I3  or  z  is  the  spanwise  direction;  u, 
V  and  w  (or  uj,  U2  and  U3)  are  the  velocity  components  in  the  coordinate  directions]. 

The  filtered  Navier-Stokes  and  continuity  equations,  which  describe  the  evolution  of  the 
large,  energy-carrying  eddies,  can  be  obtained  by  applying  the  filtering  operation  to  the 
incompressible  Navier-Stokes  and  continuity  equations  to  yield 
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in  which  a  reference  length  and  velocity  scale  are  used  to  make  Ui,  p,  Xi  and  t  dimensionless. 
The  effect  of  the  subgrid  scales  appears  in  the  residual  stress,  Tjj  =  uiuj  —  UiUj,  which,  in 
the  present  work,  was  parameterized  by  the  modified  Smagorinsky  model  [8]: 


Tij  =  -I-  6,jq]jZ  (4) 

where  6ij  is  Kronecker  delta,  and  =  Tkk  is  the  sub^rid  scale  energy  (which  is  added  to 
the  pressure),  Sij  —  {dui/dxj  +  duj/dxi)/2  is  the  large-scale  strain-rate  tensor,  and  v-r  is  an 
eddy  viscosity  given  by 


UT  =  [iCsiAxAyAzf^  (1  -  y/2S^j,  (5) 

where  a  superscript  indicates  a  quantity  made  dimensionless  by  the  kinematic  viscosity 
1/  and  the  shear  velocity  =  {t^I where  is  the  wall  shear  stress  and  p  the  fluid 
density.  The  scaling  factor  7  =  {Hi  —  H)/{Hi  —  Ht)  (in  which  H  =  5*/d  is  the  shape  factor, 
6*  is  the  displacement  thickness,  9  the  momentum  thickness,  and  the  subscripts  I  and  t  refer 
respectively  to  laminar  and  fully  developed  turbulent  flow)  was  introduced  to  decrease  the 
dissipation  by  the  subgrid  scales  during  the  early  stages  of  transition  [8]. 

The  scaled  Smagorinsky  model  described  above  was  applied  to  the  simulation  of  transition 
in  a  plane  channel  at  Re  =  8,000  (based  on  the  channel  half-width  8  and  on  the  centerline 
velocity  in  laminar  flow,  Uc).  Initial  conditions  consisted  of  the  parabolic  mean  flow,  on 
which  a  2D  Tollmien-Schlichting  (TS)  mode  of  2%  amplitude  and  a  3D  TS  mode  of  0.02% 
amplitude  were  superimposed.  The  initial  conditions  and  Reynolds  number  matched  those 
of  the  direct  simulation  described  in  [12].  The  governing  equations  (2-3)  were  integrated 
in  time  using  a  Fourier-Chebyshev  collocation  method  [13].  The  sharp  Fourier  cutoff  filter 
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was  employed  [8].  The  model  constant  (7,  was  set  equal  to  0.1,  and  the  values  5/2  and  1.7 
were  used  for  Hi  and  Ht-  The  numerical  results,  at  least  in  the  fully  developed  turbulent 
regime,  are  not  expected  to  depend  much  on  the  values  of  Hi  and  Ht-  turbulent  statistics 
are  insensitive  to  changes  in  C,  of  the  order  of  20%  (2].  A  mesh  using  up  to  48x64x48  grid 
points  was  used.  The  large-eddy  simulation  results  will  be  compared  with  those  of  a  coarse 
DNS  which  used  the  same  number  of  points  as  the  LES,  but  no  model,  and  with  the  finely 
resolved  DNS  [12],  which  used  216x216x192  grid  points.  Although  the  LES  and  the  DNS 
codes  were  essentially  the  same,  and  thus  the  CPU  time  per  point  and  per  step  was  equal, 
the  smaller  number  of  points  used  in  the  LES  resulted  in  savings,  in  terms  of  CPU  time 
required  to  calculate  the  full  transition,  of  over  90%. 

3.  Results  and  discussion 

The  time-development  of  the  wall  shear  is  shown  in  Figure  la  (in  the  following,  time  is 
normalized  by  SjUc,  velocities  by  Uc,  and  lengths  by  f;  moreover,  <  .  >  denotes  plane 
averaging,  Ui  =<  Ui  >,  and  the  resolved  fluctuations  are  defined  as  u"  =  Ui  —  Ui).  During  the 
linear  stages  of  the  growth  of  the  perturbation  (t  <  160)  all  simulations  agree  very  well.  The 
two  coarse  simulations,  however,  predict  an  earlier  occurrence  of  the  second  spike  stage  than 
the  direct  simulation.  Little  difference  can  be  observed  between  the  coarse  simulations  for 
t  <  180,  since  up  to  that  time  the  mean  velocity  has  been  affected  little  by  the  perturbation. 
It  should  be  remarked  here  that  up  to  t  =  165,  16x16x48  grid  points  or  fewer  are  necessary 
to  resolve  the  velocity  field.  For  t  >  165,  however,  the  mesh  is  refined  in  stages,  as  increased 
resolution  is  required,  to  its  final  value.  During  laminar  breakdown,  moreover,  very  large 
velocity  fluctuations  are  observed;  to  satisfy  the  CFL  condition,  therefore,  the  timestep 
must  be  significantly  decreased.  As  a  result  of  these  two  constraints,  the  first  165  time  units 
required  only  a  small  fraction  of  the  CPU  time  required  for  the  entire  simulation:  on  a  Cray 
YMP,  3600  CPU  seconds  were  used  to  integrate  the  Navier-Stokes  equations  from  t  =  0  to 
t  =  165,  while  26,500  CPU  seconds  were  expended  to  advance  the  solution  from  t  =  165  to 
t  =  250. 

Although  up  to  t  =  175  the  perturbation  energy  level  is  small  enough  that  the  interaction 
of  the  perturbations  with  the  mean  flow  is  negligible,  so  that  the  shape  factor  (Figure  lb)  is 
close  to  its  laminar  value  and  the  eddy  viscosity  i/t  is  essentially  zero  (Figure  2),  after  that 
time  they  start  affecting  the  mean  flow,  as  indicated  by  the  rapid  rise  in  the  wall  shear  and 
in  the  decrease  in  the  shape  factor.  As  the  shape  factor  is  altered,  the  eddy  viscosity  begins 
to  increase;  at  t  =  220  it  is  approximately  40%  larger  than  in  fully-developed  turbulent 
flow.  This  is  due  to  the  large  velocity  fluctuations  that  occur  during  laminar  breakdown 
(see  below),  which  cause  oscillations  of  the  large-scale  strain  rate  tensor  of  greater  amplitude 
than  in  fully-developed  turbulent  flow.  The  discrepancy  observed  around  t  =  170  between 
coarse  and  fine  simulations  appears  to  be  caused  mostly  by  insufficient  energy  transfer  to 
the  small  scales  and  viscous  dissipation.  Although  in  this  flow  the  scaling  factor  7  causes 
underprediction  of  the  energy  transfer  from  large  to  small  scales,  for  the  flat-plate  boundary 
layer  the  opposite  effect  was  observed  [8],  perhaps  because  of  deficiencies  in  the  form  of  7 
chosen  in  the  present  work. 

As  the  flow  undergoes  laminar  breakdown,  more  small  scales  are  generated,  and  the 
dissipation  provided  by  the  increased  eddy  viscosity  becomes  significant,  and  the  LES  results 
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become  more  accurate.  The  coarse  DNS,  on  the  other  hand,  further  loses  accuracy  because 
the  viscous  stresses  cannot  provide  sufficient  dissipation.  In  particular,  the  overshoot  in 
the  wall  shear  is  captured  accurately  by  the  LES,  the  shape  factor  development  is  well 
predicted,  and  the  mean  velocity  profiles  agree  well  with  the  fine  DNS  results  (Figure  3). 
The  peak  wall  shear  computed  with  the  LES  is  within  3%  of  the  DNS  prediction,  while  that 
calculated  with  the  coarse  DNS  is  in  error  by  approximately  12%.  While  at  i  =  176  the 
Reynolds  stresses  are  overpredicted  by  the  LES  because  of  the  accelerated  development  of 
the  disturbance  (Figures  4  and  5),  during  laminar  breakdown  {t  >  185)  better  agreement 
between  LES  and  DNS  is  observed;  in  those  cases  also  the  LES  results  are  more  accurate  than 
those  obtained  from  the  coarse  DNS.  Time-averaged  velocity  and  Reynolds  stress  profiles  in 
fully- developed  turbulent  flow  (obtained  from  the  large-eddy  simulation)  are  also  shown  for 
comparison.  During  the  late  stages  of  transition  the  Reynolds  stresses  can  be  several  times 
larger  than  their  counterparts  in  turbulent  flow,  perhaps  due  to  the  highly  intermittent 
character  of  the  late  stages  of  the  transition  process,  which  results  in  large  velocity  gradients 
and  increased  turbulent  kinetic  energy  production.  The  evolution  of  the  energy  content  of 
selected  harmonics,  shown  in  Figure  6,  exhibits  similar  trends  to  those  described  above:  for 
t  >  165  the  presence  of  the  subgrid  scale  model  gives  improved  results.  The  large-eddy 
simulation  was  carried  out  into  the  turbulent  regime,  and  its  results  were  found  to  compare 
well  experimental  and  numerical  results  [14,  15]. 

4.  Conclusions 

The  present  study  indicates  that,  although  a  coarse  direct  simulation  is  capable  of  repro¬ 
ducing  the  entire  transition  process  from  laminar  to  turbulent  flow,  during  the  late  stages 
of  laminar-turbulent  transition  the  presence  of  a  residual  stress  model  is  beneficial,  and 
yields  more  accurate  evolution  of  the  perturbations,  as  evidenced  by  the  comparison  of  var¬ 
ious  quantities  of  engineering  interest  such  as  mean  velocity  profiles,  Reynolds  stresses,  wall 
shear  and  shape  factor.  The  evolution  of  single  harmonics  is  also  predicted  more  accurately 
when  a  residual  stress  model  is  used  than  when  none  is  employed. 

The  predictions  of  the  scaled  Smagorinsky  model  [8]  are  not  very  accurate  at  the  early 
nonlinear  stages  of  transition,  but  the  model  becomes  more  accurate  as  small  scales  evolve 
towards  local  equilibrium  (in  which  they  receive  energy  from  the  large  scales  and  transfer  it 
entirely  to  smaller  scales  until  viscous  dissipation  becomes  significant).  Since  the  Smagorin¬ 
sky  model  is  based  on  the  local  equilibrium  assumption  for  the  small  scales,  it  is  not  surprising 
that  the  model  fares  better  at  the  late  stages  of  transition  than  at  the  early  nonlinear  ones. 
Presumably,  for  engineering  flows  in  which  a  high  level  of  disturbances  is  present,  this  type  of 
model  could  give  satisfactory  results.  It  appears,  however,  desirable  to  develop  better  ways 
to  relate  the  length  and  velocity  scales  present  in  (5)  to  the  energy  content  of  the  unresolved 
scales. 

The  fact  that  in  the  boundary  layer  case  the  scaling  factor  7  produced  a  slowed  down 
perturbation  growth,  while  in  this  flow  accelerated  growth  is  observed,  indicates  that  the 
present  formulation  cannot  fully  account  for  the  very  complex  flow  physics.  More  accurate 
scaling  or  intermittency  factors  must  be  developed,  which  depend  strongly  on  the  local  state 
of  the  flow.  This  local  character  would  facilitate  large-eddy  simulation  of  natural  transition, 
which  usually  takes  place  through  the  development  of  localized  regions  of  turbulence  (spots). 
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in  which  the  velocity  statistics  resemble  turbulent  ones,  surrounded  by  relatively  quiescent 
areas.  A  model  which  exhibits  this  behavior  is  the  RNG  model  used  in  [9];  this  model  is 
still,  however,  based  on  local  equilibrium  assumptions;  its  results,  moreover,  were  found  to 
be  very  grid-dependent,  indicating  that  a  more  accurate  length  scale  than  that  employed 
in  [9]  is  necessary.  Alternatively,  use  of  a  transport  equation  for  the  velocity  scale  valid 
during  transition  could  also  result  in  more  accurate  prediction  of  the  residual  stresses  during 
transition.  These  approaches  are  presently  under  investigation. 


5 


References 

[1]  J.W.  DeardorfF.  J.  Fluid  Mech.  41,  453  (1970). 

[2]  P.  Moin  and  J.  Kim.  J.  Fluid  Mech.  118,  341  (1982). 

[3]  U.  Piomelli,  J.H.  Ferziger,  P.  Moin  and  J.  Kim.  Phys.  Fluids  A  1(6),  (1989). 

[4]  L.  Schmitt,  K.  Richter  and  R.  Friedrich.  Finite  Approx,  in  Fluid  Mech.  (E.H.  Hirschel, 
ed.).  Vieweg,  232  (1985). 

[5]  U.  Piomelli,  P.  Moin  and  J.H.  Ferziger.  AIAA  Paper  No.  89-0375,  (1989). 

[6]  K.  Horiuti.  J.  Phys.  Soc.  Japan  55,  No.  5,  1528  (1986). 

[7]  K.  Dang  and  V.  Deschamps.  In  Numerical  Methods  in  Laminar  and  Turbulent  Flows, 
(C.  Taylor,  W.G.  Habashi  and  M.M.  Hafez,  eds.).  Pine  Ridge,  p.  423  (1987). 

[8]  U.  Piomelli,  T.A.  Zang,  C.G.  Speziale  and  M.Y.  Hussaini.  Phys.  Fluids  A  2(2),  257 
(1990). 

[9]  U.  Piomelli,  T.A.  Zang,  C.G.  Speziale  and  T.S.  Lund.  In  Instability  and  Transition, 
(M.Y.  Hussaini  and  R.G.  Voigt,  eds.),  vol.  II,  Springer-Verlag,  p.  480  (1990). 

[10]  J.  Smagorinsky.  Monthly  Weather  Review  91,  99  (1963). 

[11]  V.  Yakhot  and  S.A.  Orszag.  J.  Sci.  Computing,  1(1),  3  (1986). 

[12]  T.A.  Zang,  N.  Gilbert  and  L.  Kleiser.  In  Instability  and  Transition,  (M.Y.  Hussaini  and 
R.G.  Voigt,  eds.),  vol.  II,  Springer-Verlag,  p.  283  (1990). 

[13]  T.A.  Zang  and  M.Y.  Hussaini.  In  Nonlinear  wave  interactions  in  fluids,  (R.W.  Miksad, 
T.R.  Akylas  and  T.  Herbert,  eds.),  ASME,  p.  131  (1987). 

[14]  R.B.  Dean.  Jour.  Fluids  Engng.  100,  215  (1979). 

[15]  J.  Kim,  P.  Moin  and  R.D.  Moser.  J.  Fluid  Mech.  177,  133  (1987). 


6 


16 


>- 

"D 

Z) 

"D 


Figure  1:  Time  evolution  of  the  wall  shear  and  shape  factor.  .  Coarse  dir 

- large-eddy  simulation;  A  fine  direct  simulation  [12].  (a)  Wall  shear;  (b] 
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Figure  3:  Plane- averaged  mean  velocity.  .  Coarse  direct  simulation;  - 

simulation;  A  fine  direct  simulation  [12],  □  fully- developed  turbulent  flow, 
(b)  t  =  200;  (c)  t  =  220. 


-  large-eddy 
a)  t  =  176; 
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Figure  5:  Plane- averaged  shear  Reynolds  stress  <  u"v"  >.  .  Coarse  direct  simulation; 

- large-eddy  simulation;  A  fine  direct  simulation  [12],  □  fully- developed  turbulert  flow. 

(a)  t  =  176;  (b)  t  =  200;  (c)  t  =  220. 
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Figure  6:  Evolution  of  selected  harmonics.  .  Coarse  direct  simulation; - large-eddy 

simulation;  A  fine  direct  simulation  [12].  (a)  (0,1)  mode;  (b)  (1,0)  mode;  (c)  (1,1)  mode;  (d) 
(2,2)  mode. 
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